Ground state of a confined Yukawa plasma 
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The ground state of an externally confined one-component Yukawa plasma is derived analytically. 
In particular, the radial density profile is computed. The results agree very well with computer 
simulations on three-dimensional spherical Coulomb crystals. We conclude in presenting an exact 
equation for the density distribution for a confinement potential of arbitrary geometry. 
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I. INTRODUCTION 

Plasmas in external trapping potentials have been at- 
tracting increasing interest over the last few years in 
many fields, including trapped ions, e.g. 0, 0], dusty 
plasmas, e.g. 0- El and electrons and positrons in 
Penning traps, see e.g. Q for an overview. Among the 
main reasons is that, in these systems, it is relatively 
easy to realize strong correlation effects in charged par- 
ticle systems. Probably the most spectacular manifes- 
tation of these effects is Coulomb liquid behavior and 
crystal formation which has been found in various ge- 
ometries. In particular, the ion crystals and the recently 
observed spherical dust crystals or "Coulomb balls" 
have triggered intensive new experimental and theoreti- 
cal work, e.g. The shell structure of these crys- 
tals, including details of the shell radii and the particle 
distribution over the shells has been very well explained 
theoretically by a simple model involving an isotropic 
Yukawa-type pair repulsion and an harmonic external 
confinement potential |lOj| . 

Still, it remains an open question, how the average 
particle distribution inside the trap looks like, if it is the 
same as in the case of Coulomb interaction. It is well 
known that in a parabolic potential, particles interacting 
via the Coulomb potential establish a radially constant 
density profile. Here, we extend this analysis to a plasma 
with Yukawa interaction by solving a variational prob- 
lem for the ground state density (Sec. II). Then, in Sec. 
Ill we demonstrate that screening has a dramatic effect 
on the density profile giving rise to a parabolic decrease 
away from the trap center. There we demonstrate that 
the result for the density profile can be directly general- 
ized to any anisotropic confinement potential. While our 
analysis is based on a continuous plasma model on the 
mean- field level, we find (Sec. IV), by comparison with 
molecular dynamics simulations, that the results apply 
also to spherical crystals with a shell structure. 



II. GROUND STATE OF A CONFINED 
PLASMA 

We consider a finite one-component plasma (OCP) 
containing N identical particles with mass m and charge 
Q in an external potential $ with pair interaction poten- 
tial V described by the hamiltonian 
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The classical ground state energy follows from Eq. for 
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vanishing momenta and can be written as 0, lll| 

E[n] = J d 3 ru(r), 

with the potential energy density 
N - 1 



u(r) = n(r)< $(r) + 
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being a functional of the density profile n(r), and we ne- 
glected correlation contributions. The ground state cor- 
responds to the minimum of the energy J2J) with respect 
to the density profile with the restrictions that the den- 
sity is non-negative everywhere and reproduces the total 
particle number, i.e. 



d 3 rn(r) = N. 

This gives rise to the variational problem 
SE[n, /j] 



= 
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where 



E[n, n] = E[n] + n< N - / d 3 r n(r) 



(4) 
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and we introduced a Lagrange multiplier \x (the chemical 
potential) to fulfil condition The variation leads to 
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which holds at any point where the density is non-zero. 
Also, Eq. (J7J is equivalent to vanishing of the total force 
on the particles separately at any space point r, cf. 

sec. nun 

Equation Q is completely general, applying to any 
pair interaction V and confinement potentials of ar- 
bitrary form and symmetry, see Sec. IIII CI Of par- 
ticular interest is the case of an isotropic confinement 
$(r) = <&(r), which leads to an isotropic density distri- 
bution n(r) = n(r) = n(r)0(i£— r) the outer radius R of 
which is being fixed by the normalization condition ifljl 
which now becomes J Q drr 2 n(r) = N/Att. 



III. DENSITY PROFILE OF A YUKAWA OCP 

We now consider the case of an isotropic Yukawa pair 

potential, V(r) — ^-e' Kr which trivially includes the 
Coulomb case in the limit k — > 0. Carrying out the angle 
integration in the interaction energy in Eq. Q we obtain 
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This equation is the desired connection between the 
ground state density n(r) of the Yukawa plasma with 
the external confinement <i>(r). This integral equation 
can be solved for the density by differentiating two times 
with respect to r ^| with the result 
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±Q 2 n(r) = + $"(r)-K 2 $ + K V (9) 



The yet unknown Lagrange multiplier can be obtained 
by inserting this explicit solution into Eq. JSJ, which is 
then treated as an equation for [i, with the result 



R<f>'{R) 
1 + kR' 
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Finally, the outer radius R limiting the density profile is 
calculated from the normalization with the result 

O 2 O 2 
- 15— (N - 1) - 15— k(N - 1)R 
a a 

+ 15i? 3 + 15kR 4 + 6k 2 R 5 + k 3 R 6 = 0. (13) 

This equation has four complex and two real solutions, 
only one of which is non-negative and thus constitutes the 
unique proper result entering Eq. I|12(l . In the Coulomb 
limit, Eq. (|L L|) reduces to the familiar result of a step 
profile, 
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where the outer radius is given by 
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which is fixed by the number of particles and the con- 
stant density, the latter being controlled by the curva- 
ture a of the potential. In the right part of Eq. i|15|) we 
introduced the length scale ro = ^/2Q 2 /a, which is the 
stable distance of two charged particles in the absence of 
screening and which will be used below as the proper 
unit for lengths, screening parameter and density. Note 
that Eq. (|15|l holds also for a weakly screened Yukawa 
plasma with kR <C 1. 



A. Parabolic confinement potential 




For the frequently encountered case of a parabolic ex- 
ternal potential, $(r) = ^r 2 , we obtain for the density 
from Eq. © 

n(r) 
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where the constant c is given by 

R 2 k 2 3 + kR 
2 1 + kR' 



c = 3 



(11) 
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FIG. 1: Plasma cloud radius R [positive real solution of 
Eq. I13H 1 for a parabolic confinement potential $(r) = j r 2 
and Yukawa interaction with screening parameter k. 

In the other limiting case, kR ^ 1, the radius has 
the asymptotics kR ps (nr ) 3 (N — l)] 1 / 5 — 1. In gen- 
eral Eq. I|13l) cannot be solved for R explicitly. However, 
a general analytical result can be found by noting that 
all parameters entering Eq. l|13f) combine into only two 
parameters, x = (rero) 3 (./V — 1) and y = kR. Introduc- 
ing these paramters into Eq. i|13[l . an explicit solution is 
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found for the inverse function x(y), which can be written 
as 



x(y) 



2y A y A + 6y 2 + 15y + 15 
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Fig.^shows the result for the dimensionless radius kR = 
y of the plasma cloud, i.e. the solution of Eq. l|13fl which 
for all values of k and N is given by a single curve. 

With the result for R the constant c, which is propor- 
tional to central density, can be computed from Eq. (|12|l , 
and the complete density profile, Eq. Qllfl. is found. The 
results are shown in Fig. [21 for four particle numbers be- 
tween TV = 100 and N = 2000. One clearly recognizes 
the inverted parabola which terminates in a finite density 
value, i.e. in a discontinuity, at r = R. With increas- 
ing N, the density increases continuously at every space 
point and, at the same time, extends to higher values R. 
Thereby the density profile retains its shape. 
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FIG. 3: (Color online) Dependence of the central density n(0) 
(red dotted line), density at the boundary, n(R), (blue dashed 
line) and compression n(0)/n(R) of the plasma (green full 
line) as a function of particle number and screening parame- 
ter. 




FIG. 2: (Color online) Radial density profile for a parabolic 
confinement potential <3>(r) = ar 2 /2 and a constant screening 
parameter nro = 1 and four different particle numbers N 
shown in the figure. For comparison, also the result for ftro = 
0.3 and TV = 2000 is shown by the dashed line. 

On the other hand, when the plasma screening is in- 
creased, at constant N, the density profile changes dra- 
matically, compare the two curves for N = 2000. Increase 
of k leads to compression of the plasma: the radius R de- 
creases, and the absolute value of the density increases, 
most significantly in the center. This compressional be- 
havior is shown in Fig. |3 cf. the full green line showing 
the ratio of the inner to outer densities of the plasma. 

The dependence on k is analyzed more in detail in 
Fig. 0] below for a fixed particle number N = 2000. In 
the case of Coulomb interaction, k = 0, we recover the 
constant density profile (|14fl . On the other hand, in the 
case of a screened potential, the density decays parabol- 
ically with increasing distance from the trap center, cf. 
Eq. i|ll|) . Also, the density discontinuity at r = R is soft- 
ened compared to the Coulomb case, and the step height 
increases. 
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FIG. 4: (Color online) Radial density profile, solutions of 
Eq. lllll . of a three-dimensional plasma of N = 2000 for 
four values of the screening parameter (lines), from bottom 
to top: k — (red), kxo = 0.3 (yellow), kxq = 1 (green), 
Kro = 3 (blue). Symbols denote molecular dynamics results 
of a plasma crystal for the same parameters where the aver- 
age density at the positions of the shells is shown, for details 
see discussion in Sec. lIVBl 



B. Force equilibrium 

Besides minimizing the total energy, cf. Eq. Q, the 
ground state density profile can be obtained from the 
condition of a local equilibrium of the total force (inter- 
nal plus external ones) at each point where the density 
is non-zero. On the other hand, the shape of the radial 
density profile and its change with variation of n is di- 
rectly related to a change of the force equilibrium. Here, 
we demonstrate this for the special case of a spherically 
symmetric confinement potential $(r). The forces can 
be directly derived from Eq. by taking the gradient 
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in radial direction 



&(r)=F < (r) + F > (r), 



(17) 



which means that, for any spherical layer at a distance r 
from the center, the external force F$(r) — — <f>'(r) which 
acts towards the center is balanced by the internal force 
due to the Yukawa repulsion between the particles. The 
internal force contains two parts where -F< arises from 
the action of all particles inside the given layer, r' < r, 
and acts outward, whereas F > results from the action of 
all particles located outside r' > r, and acts inward, 



F<(r) = 4^^Q 2 - 



N- 1 „,1 



1-1 I x 

r \ nr 



F>(r) = A-K——Q 1 - - cosher) + 

N r \ nr 



dr' r n(r ) sinh(Kr'), 



sinh(«x) 



dr'r'h{r')e Kr 



This force balance can be used to obtain the ground state 
density profile. Alternatively, we can use the computed 
profile to analyze the two internal force contributions and 
their dependence on k. 

Consider first the limit of weak screening, kR <C 1. 
Then the forces approach the Coulomb case and, in the 
case of a constant density profile 114|) . 



N - 1 O 2 



N r 2 



(18a) 
(18b) 



with N < = n c 47rr 3 /3 being the particle number in the 
inner region. This means, the force is repulsive and in- 
creases linearly with r and exactly compensates the linear 
external force F$(r) = —ar for all values r < R. 

In the general case of finite screening the outer force 
F > (r) does not vanish, cf. Eq. (|18a|) . Since its direction 
is always towards the center, the force F<(r) has to in- 
crease simultaneously, in order to compensate the com- 
bined effect of F$ (r) and F> (r) . This effect increases 
continuously with increasing k which is directly verified 
by evaluating the expressions in Eq. (|18afl . 



C. Generalization to arbitrary confinement 
geometry 

The result for the density profile in an isotropic con- 
finement, Eq. can be easily extended to arbitrary 
geometry. For this purpose we use the text book result 
that the charge density corresponding to the Yukawa po- 
tential is QS(r) — QK 2 e~ KT jr. This allows us to rewrite 



the Poisson equation as 



(A — k 2 



showing that the Yukawa potential is the Green's func- 
tion of Eq. (fB?f . This fact can be used in Eq. Q for the 
case of a confinement potential <!> of arbitrary geometry 



$(r)-/x=-l— ^ Q 2 I dVn(r') e '" 



N - 1 



(20) 



|r — r'| 

to get the explicit result for the density profile 

Q 2 n(r) = A$(r) - K 2 $(r) + «V (21) 



IV. DENSITY PROFILE OF CONFINED 
COULOMB AND YUKAWA CRYSTALS 

So far we have considered the model of a continuous 
density distribution n(r). On the other hand, the ground 
state of a confined spherically symmetric system of dis- 
crete point-like charged particles is known to have a shell 
structure as was demonstrated for dusty plasmas in Ref. 
0- It is, therefore, of interest to verify if such a shell 
structure can be derived from our starting equation @ 
for the total energy and to compare our results to the 
radial density distribution in such Coulomb or Yukawa 
balls. 



A. Derivation of a shell model for a trapped finite 
Yukawa plasma 

The concentric shells observed in spherical trapped 
Coulomb crystals have led to the proposal of simple an- 
alytical models, cf. e.g. 0,0,0,0]. Such a model 
for a trapped one-component plasma is trivially derived 
from the total energy expression © by inserting for the 
density the ansatz 



n s (r) = n s (r) 



L 



^ AttRI 



S(r-R v ), 



(22) 



which describes L concentric shells of zero thickness with 
N v particles on shell v with radius R v and N v +£ 
Y . where £ denotes the number of particles in the trap 
center (zero or one) 0, H^ |. As a result, we obtain for 
the total ground state energy of a Yukawa plasma in an 
isotropic general confinement potential $ 



E s (N;k) =J2N V \$(R„) +Q 



R v 



wMkRv) N v -1 
kR~„ 2 



sinh(K_R A1 ) 

K-Rn, 



N, 



This is essentially the Yukawa shell model of Totsuji et 
al., Ref. 0| where, however, the finite size correction 
factor (N v — 1)/N V in the intrashell contribution and the 
term £ are missing. In the Coulomb limit, k — > 0, the 



(19) 

result simplifies with e 



kR v 



1 and 



sin kRv 

KRv 



1, and 
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we immediately recover the Coulomb shell model of Hasse 
and Avilov [14( (plus the additional correction factor). 

A further improvement is possible by including in- 
trashell correlations 0] . The simplest model is obtained 
by replacing N v — 1 — > N v — e(N)y/N v , where e is a fit 
parameter close to one which allows to achieve excellent 
agreement with the exact ground state [TT| . An alter- 
native way to include correlations was proposed by Ref. 

B. Comparison with simulation results for finite 
Yukawa crystals 

In order to compare the density profile n(r) of our 
continuous model with the density of discrete spherical 
Yukawa crystals, we performed molecular dynamics sim- 
ulations of the ground state of a large number of Coulomb 
balls, for details, see refs. |[l,[l{J. As an example, the nu- 
merical results for a Coulomb ball with N = 2000 which 
is large enough to exhibit macroscopic behavior [Iti Il7j 
are included in Fig. The symbols denote the average 
particle density around each of the shells. The averag- 
ing was accomplished by substituting each particle by a 
small but finite sphere, so that a smooth radial density 
profile was obtained. 

With increasing k the crosses move towards the cen- 
ter confirming the compression of the Coulomb balls ob- 
served before ^(j- Obviously, the simulation results are 
very well reproduced by the analytical density profile 
of a continuous plasma. But there are also small discrep- 
ancies in the central part which grow slowly with k. One 
reason is that, for large k, the width of the inner shells 
increases rapidly, making the comparison difficult. An- 
other possible reason could be the effect of correlations. 



V. SUMMARY AND DISCUSSION 

In summary, we have presented a theoretical analysis 
of the ground state density profile of spatially confined 
one-component plasmas in dependence on the form of the 
pair interaction. An explicit result for the density profile 
for an arbitrary confinement potential has been derived. 
In particular, for an isotropic confinement, we have found 
that screening of the Coulomb interaction substantially 
modifies the radial density distribution. In contrast to 
a bare Coulomb interaction for which the density inside 
a parabolic external potential is constant, for a screened 
interaction, a quadratic decay away form the center is 
found. 

Interestingly, while our results were derived for a con- 
tinuous density distribution (a macroscopic system) and 
with neglect of binary correlations, our analytical results 
agree very well also with first-principle simulation re- 
sults for strongly correlated Coulomb and Yukawa clus- 
ters containing several thousands of particles for screen- 
ing paramters ktq < 1. This agreement is by no means 
trivial and deserves further analysis. Further, it is very 
interesting to investigate the reason for the deviations at 
large values of the screening parameter, and analyze the 
effect of binary correlations on the density distributions 
is. These questions will be subject of forthcoming work. 
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